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ABSTRACT 

We present a time-dependent approach to the one-zone hadronic model in the case where the photon spectrum is produced by ultra- 
relativistic protons interacting with soft photons that are produced from protons and low magnetic fields. Assuming that protons are 
injected at a certain rate in a homogeneous spherical volume containing a magnetic field, the evolution of the system can be described 
by five coupled kinetic equations, for protons, electrons, photons, neutrons, and neutrinos. Photopair and photopion interactions are 
modelled using the results of Monte-Carlo simulations and, in particular from the SOPHIA code for the latter. The coupling of energy 
losses and injection introduces a self-consistency in our approach and allows the study of the comparative relevancy of processes at 
various conditions, the efficiency of the conversion of proton luminosity to radiation, the resulting neutrino spectra, and the effects of 
time variability on proton injection, among other topics. We present some characteristic examples of the temporal behaviour of the 
system and show that this can be very different from the one exhibited by leptonic models. Furthermore, we argue that, contrary to 
the wide-held belief, there are parameter regimes where the hadronic models can become quite efficient. However, to keep the free 
parameters at a minimum and facilitate an in-depth study of the system, we have only concentrated on the case where protons are 
injected; i.e., we did not consider the effects of a co-accelerated leptonic component. 

Key words. Radiation mechanisms: non thermal; Radiative transfer; Galaxies: active 
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case by Protheroe & KazanasI (1 19831) and were later applied to 
jets by Si ewert et al.1 (12004 and lRevnoso et al.1 (1201 ll) . For a r 



BL Lac objects and flat-spectrum radio quasars (commonly 
termed blazars) display a spectral energy distribution (SED) that 
extends from radio to y-rays, in the shape of two broad humps - 
one in the radio to UV or X-ray frequency range and the other in 
the X-ray to y-ray range. Both components are non-thermal and 
very likely originate in the blazar's jet, which is closely aligned 
with the line of sight of the observer. While the lower-energy 
hump is widely acknowledged to result from synchrotron radia- 
tion from relativistic electrons in the emitting region, the origin 
of the higher-energy hump is subject to differing interpretations. 

According to the leptonic models, the high-energy compo- 
nent emerges from inverse Compton scattering of electrons on 
soft target photons. Those photons may b e the product of syn- 
chrotron radiation of the s ame electro ns (Maraschi et al.ll 1992c 
Bloo m & Marscherl 119961; llnoue & Takaharal 1 1996b or an in- 
digenous populati on of external photons (lDermeretal.lll992t 
Siko raet al.l tl994t. Hadronic models, on the other hand, as- 
sume that protons, which are accelerated along with elec- 
trons, contribute much of the high-ene r gy component via 
photopion interactio ns (Mannheiml Il998t IMiicke et all 120031: 
Bottc her et al.l l2009h . The target photons may either be pro- 
duced inside the jet via synchrotro n emission from a co- 
accelerated populatio n of electrons (LM annhei ml 1 19931). orig- 
inate outside the je t (Protheroe 1997; Bednarek & Protheroe 
1 1999t lAtovan & Dermerl 1200 ll) . or be produced via syn- 
chrotron emission from the pr otons themselves (lAharonianl 
2000; iMucke & Pro theroe 2001aj). In addition proton-proton in- 
teractions were taken into account for the first time in the AGN 



cen t review of the f eatures of leptonic versus hadronic models 
see iBoettcheil d2012l) . Although similar models a r e used to study 
gamma ray bursts dBottcher & Dermerl 1 1998t iKazanas et all 
1 20021; lMa stichiadi s~& KazanasI l2006t lAsano & Inoud 120071; 
ksanoetalJ l2009t iMastichiadis & KazanasI l2Q09h a nd y- 
ray emitting compact binary systems ([Romero et al. 2003; 



Pare des et al I 120051 iRomero et al.l 120051) . the former demand 
large magnetic fields, while the latter demand high proton den- 
sities, making inelastic pp-collisions more important than pho- 
tohadronic interactions. However, it should be noted that both 
photo-pair and photo-meson production are included in recent 
models of compact binaries, and they are shown to be signif- 
icant when proton a nd electron luminosities are near equiparti- 
tion dVila et al.l2012l) . The processes and range of variables used 
in this work are, instead, tailored to photon-dominated systems 
such as blazars, although they can indice similar effects in other 
systems, and thus are qualitatively relevant to them. 

While leptonic models have often made use of time- 
dependen t codes to solve the kinetic equations of electrons and 



photon s (Mastichiadis & Kirk] [T997I; iKrawczynski et al.l l2002t 
iBottcher & Chiang! 120021: iKatarzvris ki et al.l 120051) . hadronic 
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models have been very difficult to investigate in a similar man- 
ner, owing to the greater complexity of the particle interac- 
tions involved and the resulting time-consuming nature of sim- 
ulations. First efforts dMastichiad is & Krrklll995|) used simple 
<5-functions for the secondary production, while the effects 
of Bethe-Heitler pair production were c alcula ted in a time- 
dependent manner in Mastichiad is et al.l (120051) who used the 
Monte Carlo results of Protheroe & Johnson (1996) to model 
the process. In the present paper we extend this method to in- 
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corporate photopion interactions in detail. For this we have re- 
quired the use of the SOPHIA Mon te-Car lo event generator for 
hadronic interactions (Muc ke et al.l 120001) . whose results have 
been adapted into a time-dependent code. 

The paper is structured as follows. In §2 we present the time- 
dependent kinetic equations for all five particle types and show 
the coupling that emerges when all leptonic and hadronic pro- 
cesses are taken into account. In §3 we discuss the integration of 
the SOPHIA code's results into our model and its verifications 
through energy-loss tests. In §4 we present various results from 
monoenergetic proton injection in a magnetic field, noticing the 
transition in importance from photopair to photopion when the 
injected proton energy is increased, and the transition from lin- 
ear to non-linear proton cooling when the injected proton com- 
pactness is increased; also the effects of time variability in pro- 
ton injection. We move on to power-law proton injection in §5, 
while in §6 we show some characteristic time-dependent cases. 
Finally, in §7 we conclude with a summary and make some re- 
marks on the results of the previous sections. 

2. The kinetic equation approach 

2.1. Setup of the problem 

Relativistic protons in compact sources undergo a series of pro- 
cesses that can be quite complicated to model. Thus protons in- 
teract with photons creating pions, neutrons and neutrinos (all 
from photopion interactions), and also electron-positron pairs 
(from photopair). Since the magnetic field in the cases exam- 
ined in this paper is low (B < 10 G) and the maximum pro- 
ton Lorentz factor does not exceed y p = 10 8 , the pions and 
their resulting muons decay almost instantaneously, as compared 
with any other relevant timescale of the problem, creating pho- 
tons and more electrons, positrons, and neutrinos. For higher 
magnetic fields, those pions and muons would suffer significant 
energy l osses through synchrotron radiation before dec aying 
(see e.g. iMiicke & Protherod d2001bl) . iMiicke et all d2003l) and 
iRevnoso & Romerol d2009l) ) The electrons and positrons cool by 
producing photons, which can, in turn, interact with protons, in 
turn creating more secondaries. At the same time, neutrinos will 
escape, while neutrons, depending on the opacity, can either es- 
cape, interact, or decay inside the source. Since the target ma- 
terial in the studied environment has low density, proton-proton 
interactions and bremsstrahlung can be neglected. To treat this 
complicate d system, we used the kinet ic equation approach as 
descr ibed in Mastichiadi s & Kirkld 1995b : henceforth MK95, and 
iMastichiadis et al.l d2005l) : henceforth MPK05. However, in the 
present treatment we have extended the set of equations by in- 
cluding, in addition to the three for protons, electrons/positrons, 
and photons, two more equations, one for neutrons and one for 
neutrinos. The target photons and all other particle distributions 
are again taken to be isotropic. 

As in MPK05 we can write the kinetic equations for a homo- 
geneous region containing relativistic hadrons and leptons in the 
compact form 



^ + L, + e , = o. 

at 



(1) 



Here, index i refers to protons (denoted by 'p'), elec- 
trons/positrons ('e'), photons ('y'), neutrons ('n') and neutrinos 
(V). The functions to be determined are the dimensionless dif- 
ferential number densities of the five species, which can be pro- 
duced/injected in the source through the operators Q, and de- 
stroyed/escape through the operators L,-. They are normalised 



from the ordinary differential number densities n, in the follow- 
ing way: 



niiyt, t)Aji = o-jRhiiEi, t)AE t with 



7i 



llliC 1 



(2) 



where R is the radius of the source and crj the Thomson cross 
section. The photon and neutrino energy is normalised with re- 
spect to the electron rest mass. Time is also normalised, to the 
light-crossing time of the source f cr = R/c, so t — ct/R, where i 
is the actual time. 

The processes that we have included are 

1 . Proton-photon pair production, 

2. Proton-photon pion production, 

3. Proton synchrotron radiation, 

4. Electron synchrotron radiation, 

5. Synchrotron self absorption, 

6. Electron inverse Compton scattering, 

7. Photon-photon pair production, 

8. Electron-positron annihilation, 

9. Compton scattering of photons on cool pairs, 
10. Triplet pair production. 

The modelling of t he leptonic processes ( 4-9) h as been dis- 
cussed in MK95 and [Mastichia dis & Kirk! (1 19971) . Photopair 
production (1) has been extensively discussed and modelled 
in MPK05. Modelling of proton synchrotron radiation (3) was 
done in a manner analogous to electron sync hrotron. Triplet pair 
production (10) was mod elled according to (Mastichiadis et al.l 
119861: Irvlastichiadisl 1 1 99 1 1) . However, the basic improvement of 
the present paper is for photopion production (2). This process 
was modelled in MK95 with (5-function approximations and 
without taking neutrons and neutrinos into account in the injec- 
tion and energy loss terms. Here we revisit th e process by doing a 
modelling based on the results of SOPHIA (Mu cke et al 1 2000). 
which we discuss in the next section. A simil ar modelling of the 
above processes was recently undertaken by IVievro & Romerol 
but in a strongly magnetised and "dirtier" environment, 
with a stronger emphasis on proton-proton interactions and a 
weaker emphasis on the precise modelling of phot opion interac- 
tions, which uses the approximations discussed in lRomero et al.l 

d20ioh . 

By including the various relevant terms, the kinetic equations 
for each species become 
- Protons 



at p p 



Electrons 



+ Ll 



Y n| l /-jphotopion , ^ndecay^^ 



Q? + Ql 



+ <2 



~~dt 



+ zr + m 



QT + Ge H + QV + Of otopion + 6e pp + cay ; 

- Photons 



(4) 



-7T- + — 

8t ty 



+ 1% + !%* = 



QT + Q p ; yn + + Q y nn + Q y 



photopion m 



(5) 
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- Neutrinos 



dll v n v _ ^-jphotopion ^ndecay, 
Ot f esc 



- Neutrons 



, photopion ^ ndecay _ ^-,photopion 

— — + L n + L n + - — y n 

at teac 



(6) 



(7) 



The operators are labelled according to the processes that 
produce them, i.e. proton-photon pair production (BH), proton- 
photon pion production (photopion), electron synchrotron (syn), 
proton synchrotron (psyn), synchrotron self absorption (ssa), in- 
verse Compton scattering (ics), photon-photon pair production 
(yy), triplet pair production (tpp), neutron decay (ndecay) and 
electron-positron annihilation (ann), while (ext) refers to exter- 
nal injection. f p esc and f e esc are the escape times for protons and 
electrons, respectively, while f esc is the escape time for neutrally 
charged particles. Compton scattering of photons on cool pairs 
is approximated by multiplying the photon escape term, f r e sc, by 
the factor 



fy.esc = fee X (1 + H(l - x)T T /3) 1 



(8) 



where tj is the Tho mson optical depth and H{\ - x) is the 
Heavyside function (Lightman & Zdziarski] 1 19871) . with x = 
hv/m e c 2 the dimensionless photon energy. 



their minimum energies forming a tail that can extend over sev- 
eral orders of magnitude. Although these distributions have re- 
cently been approximated analytically (see Kelner & Aha ronianl 
2008), we made use of them in their original form, representing 
their effect with a coefficient d(y, y'), where y is the resulting 
nucleon's energy, £ the target photon energy, and y' the initial 
nucleon's energy. Each interaction can be seen as removing the 
initial nucleon from its energy bin and creating the distribution 
of new ones, so that translates into a catastrophic loss term cou- 
pled with an injection term. The cross section <x N (y',£) photo P lon 
is related to the interaction time TN(y',£) photoplon for a single nu- 
cleon in a monoenergetic photon field of density 1cm 3 , so that 
o- N (y',£) photopion = l/(T N (y',£) photopion c), where N is the type of 
nucleon. 

Therefore, the contributions to the loss and injection terms 
in the proton and neutron kinetic equations from the above pro- 
cesses are (respectively) 

- Channels (a) and (b) 



QT 



T (y, = 2 ^p(r'^) photopion «'(y, f , y> P (y', 0« r (£ (14) 

Y4 



Lf^{y, t) = - J] <r p (y', £) phot ° p '° n » p (y, tynytf, t) (15) 

f 

i P "(y, = - Yj a- v (y',O photopion n v (y, t)n 7 (t, f) (16) 



2.2. Photopion interactions: Loss and source terms 

Photopion interactions produce a distribution of neutrons that is 
similar in energy to the outgoing protons, since the two basic 
channels 



(«) 
(h) 



p+y 
p + y 



p + if 



n + n 



(9) 
(10) 



are about equally probable. Neutrons can escape the emission 
region essentially in a crossing time, but they are susceptible to 
two processes along the way. They can interact with ambient 
photons in much the same way as protons, leading to a mirror 
image of the previous channels 



(c) 



id) 



n + y 



n + y ■ 



n + n 



p + 71 , 



(11) 



or they can decay into protons, with a mean lifetime of t = 
881.5 ± 1.5s (Nakamura . K. et al. (Particle Data Group)l[?2 010) 



(«) 



p + e + v e 



(13) 



- Channels (c) and (d) 



er™(y,0 = ^cr n (y',a photoplon £ /(y^,y')«„(y , ,f)" y (^fX18) 

■ ^ o- n (y', £) phot ° p '° n « n (y, 0«y(£ (19) 

f 

■2^(y',a photop,on "n(y,0» r (^0 



r yn— >rur 



(y,0 



. yn—>pn 



(y,t) = - 



(20) 



(12) G"(T' f ) = Z ^n(y'^) photoplo 'V(y, f , y>„(y', 0»y(£ t). (21) 

Y4 



Neutron decay produces complementary loss and injection 
terms for the neutrons and protons, respectively: 



Lr p (y,t) 



In channels (a) - (d), the initial nucleons lose a portion of 
their energy that depends on the inelasticity of the interaction, 
k p , but that is not a uniform process. SOPHIA simulation re- 
sults show a distribution of nucleons after each interaction, with 
their maximum energies close to that of the initial nucleon and Qv ~ 



" n (y, t) 

yr 



»n(y, t) 



(22) 



(23) 
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However, since neutrons are produced predominantly with 
energies close to that of the protons undergoing photopion inter- 
actions, their Lorentz factors are high enough to ensure that neu- 
tron decay occurs mostly outside the source. We ignore the en- 
suing interactions of the products. Neutrinos are only produced 
in photopion interactions, and they freely escape the emission 
region, with an escape time equal to the crossing time of the 
source. In the code the neutrino flux is the sum of all flavours 
since we may assume that there will be complete mixing after 
propagation. Since neutrinos escape without interaction, their 
spectra is proportional to their injection spectrum. We neglect 
neutrinos from neutron decay. 



photons - 
neutrinos 
electrons 
protons 

neutrons 



3. Modelling photopion interactions 

3.1. The SOPHIA code 

At low interaction energies, photohadronic processes are dom- 
inated by the A(1232) resonance, so their cross section can be 
approximated by a 5-function. However, such an approximation 
neglects processes that become dominant for higher energies and 
are still important in lower ones. To remedy this, the SOPHIA 
Mont e-Carlo event generator was developed (see Muck e~et alJ 
2000), taking into account the following interaction processes: 
resonance production, direct pion production and diffractive and 
non-diffractive multipion production. 

3.2. Modeling SOPHIA results 

The photopion spectra are calculated indirectly from the results 
of this SOPHIA Monte Carlo code. Protons of three distinct en- 
ergies (y p = 10 ,7p = 10 ,7 P = 10 12 , where the subscript p 
refers to protons) are allowed to interact with isotropically dis- 
tributed target photons of energies ranging from £ = 10 to 
£ = 10° 1 in logarithmic steps of 0.1. This provides us with inter- 
action rates and energy distributions of all secondary particles, 
including protons (since the original protons are treated as hav- 
ing sustained catastrophic losses, giving rise to a new distribu- 
tion over a wide range of energies) for those three specific ini- 
tial proton energies. Specifically, we obtain an energy loss term 
for protons, and injection terms for photons, electrons, protons, 
neutrons, and neutrinos. The last result from pions and muons, 
which are assumed to decay instantaneously. 

The energy grids used in our program are equally spaced 
in the logarithm of y, with a resolution of ten bins per decade 
for protons and electrons and five bins per decade for photons. 
Furthermore, the lowest grid point for both protons and electrons 
is y m in = 10 01 , while for photons it is taken to be x m \ n = by^ n , 
since the softest photons are assumed to be produced by electron 
synchrotron radiation, b = B/B c being the magnetic field in units 
of the critical field B c = m 2 e c 3 /(eK) = 4,414 X 10 13 G. 

Each particle's injection term is calculated as the energy dis- 
tribution term from SOPHIA over the interaction time for the 
initial proton and photon energies, times the number densities of 
those initial protons and photons, all normalised to our dimen- 
sionless units, 

/-.photopion , 

d{y,%,f)m{y',t)ny{t,i) 
Q (y, t) = ——r (24) 

where d(y,^,y') is the energy distribution term, as described in 
sect. 12.21 TN„{y' , £) is the interaction time, N the type of nucleon, 
i the particle type, and y is substituted by x for photons and 
E v /m e c 2 for neutrinos. 



-22 I 1 1 1 i 1 I 1 

7 7.5 8 8.5 9 9.5 

logx, logy e 

Fig. 1. Test distributions from the interaction of monoenergetic protons 
(y p = 10 6 1 ) with photons (£ = 10~ 2 ' 8 ) as produced by the numerical 
code. Units in the y axis are arbitrary. These have been checked against 
the SOPHIA results and are in very good agreement. 



Since the produced spectra depend on the product y p • for 
any proton of energy y p interacting with a photon of energy £ 
the effects will be identical to those of a proton of energy y p 
interacting with a photon of energy = £y p /y p , the only differ- 
ence being that the energies of the secondaries will be shifted by 
y p /y p . Therefore, we can use the three distinct proton energies 
(y p = 10 61 ,y p = 10 91 ,y p = 10 12 ) to calculate the spectra for 
protons of all intermediate or lower energies with this scaling 
approach. These three energies were chosen for technical rea- 
sons (in the present paper injected protons have y p <sc 10 12 ). It 
comes at an acceptably small loss of accuracy but at a gain in 
computing time and memory. 

When the initial particle is a neutron instead of a proton, 
we can use the same data but with switched labels; neutrons in- 
stead of protons, and all particles switched with their antipar- 
ticles. Comparisons with SOPHIA results from neutron-photon 
interactions show this method to be sufficiently accurate, since 
the error in energy of the switched distributions is no higher than 
3.5%. 

3.3. Energy losses 

To test the additions to the code, we performed both single- 
interaction checks and full runs for various cases of monoener- 
getic protons injected into black body radiation fields. In Fig.Q] 
we show the distributions resulting from interactions of monoen- 
ergetic protons of Lorentz factor y p = 1Q 61 with photons with 
£ = 10~ 2 8 . There, 95.6% of the proton's energy is conserved in 
the secondaries, whose contribution by particle type corresponds 
almost perfectly to the ones in the SOPHIA data. If we run the 
same test with photons of £ = 10~ 3 8 and £, = 10~ 02 , correspond- 
ing to the lowest and highest energy photons for which we can 
test the results directly, the energy conserved in the secondaries 
is 90% and 98%, respectively. 

In Fig. [2] we show the energy losses of protons in a black 
body radiation field (introduced merely for this test and not used 
elsewhere) from both photopion interactions and pair produc- 
tion. The relative con tribution of each process for each proton 
energy is as expected (IStanev et al. 20Q3). 
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-10 




Fig. 2. Average energy loss times of protons interacting with ambient Fig. 3. Production rate of secondary electrons for R = 3 x 10 16 cm, 
black body radiation of temperature 3 x 10 5 A^ with R = 10 15 cm in- 



troduced here purely for the purpose of checking energy losses. Long- 
dashed line corresponds to pair production, while full line corresponds 
to photopion interaction. Units in the y axis are dimensionless, nor- 
malised to the light crossing time. 



B = 1 G, t v = 0.4, f p>csc = t CI and y = 2.5 x 10 6 . Full line depicts 
the injection resulting from charged pion decay while long-dashed line 
is the photopair (Bethe-Heitler) injection. Short-dashed line depicts the 
pairs injected from yy absorption, which is negligible for the particular 
set of the initial parameters chosen. 



4. Monoenergetic proton injection 

4.1. Spectral signatures 

Adopting the standard picture of the one zone radiation model, 
we can assume a spherical volume of radius R and a tangled 
magnetic field of strength B. We consider only the case where 
protons are injected in the system and set the external injection 
of electrons equal to zero. The simplest way of treating proton 
injection in Eq. (© is to assume that protons are injected with a 
(5-function at some energy yo- This will faciliate some tests of 
the radiative signatures of the different hadronic processes dis- 
cussed in the previous paragraph. Thus we write 



6p nj (r P ) = G P ,u6(y P - yo) 



(25) 



where Qp.o is the proton normalisation, considered as indepen- 
dent of time. In the present treatment we find it more convenient, 
instead of defining Q p o, to equivalently define the proton com- 
pactness which is given by the relation (as in MK95). 



fp = *k J d7pi7p " l)Ql ? i7p) - 



With the <5-function assumption, we can write 

3m, It, 



Q 



p.O 



m p (yo - 1) 



(26) 



(27) 



As was shown in MPK05, the actual injected proton luminosity 
is related to t p by the relation 



4-nRm D c 3 

J-i p — V p . 



(28) 



It will also prove useful to define, in a similar manner, the 
compactnesses 



m P f 

3m e J 



dy p (y p - l)Ll(r P ). 



(29) 



where i refers to each process that creates a loss term for protons 
in Eq. 0). 



In addition to yo, l p , and tp >esc , one needs to specify initial 
conditions for the five unknowns to fully determine the sys- 
tem. Without loss of generality, we can assume that at t = 0, 
«p(7p»0) = « e (y e ,0) = n Y (x,0) = « n (7n,0) = «vCE v ,0) = 0. 
Then we can integrate the system forward in time. 

Therefore, qualitatively speaking, one expects that for t > 0, 
protons will start accumulating in the source. At the same time, 
according to Eq. © protons will lose energy by synchrotron and, 
possibly, by photopair and photopion production, while a frac- 
tion will physically escape at a rate t p l esc from the source region. 
In this case secondaries will be created making Eqs. © to © 
relevant for the evolution of the system. 

As a first example, we show the results obtained when R = 
3 x 10 16 cm, B - \ G, l p = 0.4, f p>esc = f cr , and y = 2.5 x 10 6 . 
The magnetic field strength is such as would b e expected at a dis- 
tance o f roughly 0. lpc from the base of the jet ( Komiss arov et all 
(l2007h . lO' Sullivan & Gabuzdal (l2009l) ). The parameters were 
chosen in such a way that the radiated photons will cause mini- 
mal losses on the protons, therefore the proton steady-state dis- 
tribution is very close to the one derived from Eq. © when all 
loss terms are ignored, i.e. 



3m e ^ P f P , eS c 

« P (r P ) = r<ny P - yo)- 

m p yo-l 



(30) 



This simple form of the proton distribution function faciliates 
the investigation of some important points regarding the injected 
electron and radiated photon spectrum. For this reason we do not 
consider the presence of an external photon distribution. In this 
and the next test cases, the proton synchrotron photons serve as 
targets. 

At low compactnesses, as we are considering here, electrons 
are injected into the system as secondaries mainly through two 
channels (i) from photopair and (ii) from photopion. (At higher 
compactness photon-photon pair production becomes important 
but we can neglect this for the time being.) These two processes 
are in direct competition with each other and their relative im- 
portance depends on such parameters as yo and the soft photons 
which serve as targets. Figure [3] plots the injection functions of 
these two processes for the initial parameters given above. We 
note that the two distributions have different characteristics. The 
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log x 

Fig. 4. Steady state MW spectrum of photons resulting from monoen- 
ergetic proton injection with minimal cooling. Parameters are as in the 
previous figure. Full line is the synchrotron spectrum corresponding to 
the electrons injected from photopair (dotted line) and photopion (dot- 
dashed line) - see previous figure. Long dashed line is the proton syn- 
chrotron component while short dashed line are the y-rays resulting 
from n° decay. 



injection function of photopair electrons is broader and has a 
peak at energies y e =s y p . The injection function of photopion 
electrons, on the other hand, is flatter and peaks at much higher 
energies, of the order of y e = r] ne y p , with r\ m 150. 

In the case we are considering here, the photon spectrum will 
show four distinctive features. Two of them are connected to the 
synchrotron radiation of the injected electron populations dis- 
cussed above, while the two other features are connected to pro- 
ton synchrotron and 7r°-decay. Putting all of them in ascending 
order with frequency we have 

1. Proton synchrotron radiation: Since the proton distribution 
function is a 5-function at yo, the radiated photon spectrum 
will have a peak at x — by\ (where b - B/B c , B c = 

m 2 e c 3 /(eK)). 

2. Synchrotron radiation from photopair electrons: As stated 
above, the electron injection function resulting from pho- 
topair interactions is rather broad with a peak at y e =s y p . 
Synchrotron cooling of electrons and consequent radiation 
results in a photon spectrum with peak at x =* by\. 

3 . Synchrotron radiation from photopion electrons: In complete 
analogy to the photopair, the peak of this distribution will be 
at x =i b(j] nf ,yo) 2 . 

4. y-rays from -decay: A monoenergetic proton distribution 
produces a well defined peak at x - rj nY yo, with rj nY 350. 

It is interesting to note that the ratio of where the three first 
peaks occur is fixed, i.e. (m e /m p ) : 1 : rj^ = 5xl0~ 4 : 1 : 3xl0 4 , 
which implies that the proton synchrotron peak will always be 
about eight orders of magnitude below the synchrotron peak of 
photopion electrons. Only the y— ray peak from neutral pion de- 
cay is not connected to the other three - since it does not involve 
the magnetic field, but for all astrophysically relevant magnetic 
field values, it should always be the highest energy peak. These 
features can be seen in Fig.|4]where the multi-wavelength (MW) 
spectrum corresponding to the parameters adopted above is plot- 
ted. The inverse Compton scattering component of the electrons 
is negligible. 



4.2. increasing the injected proton energy: The transition in 
importance from photopair to photopion 

While protons will always radiate through synchrotron radia- 
tion, it is a matter of the choice of the initial parameters whether 
they will undergo substantial photopair and photopion produc- 
tion. The rate of these processes reaches a maximum when the 
energy of the collision between the peak of the proton syn- 
chrotron radiation and the relativistic proton is above their re- 
spective thresholds. If we denote the proton synchrotron typical 
energy by jc psyn = ^byl, then if 



XpsynTO < 2, 



(31) 



then obviously the proton-photon collisions cannot produce co- 
pious photopairs because most of the collisions occur below 
threshold for this process. If 



2 < -XpsynTO < — (1 + ~ — )> 
m e Zm n 



(32) 



then the collisions between jc psyn and yo are above the thresh- 
old for photopair but not for photopion; note however the syn- 
chrotron radiation from the produced photopairs can produce 
photopions, therefore we expect such a component even in this 
case, albeit at a low level. Finally, for 



m n m„ 

(1 + ~ ) < -KpsynTO 

nu 2m„ 



(33) 



we expect that protons can produce both pairs and pions directly 
in collisions with proton synchrotron photons. 

The above relations help us define two characteristic proton 
energies that are functions of the magnetic field only; i.e., 



/p.pair — I , 



1/3 



and 



7pj 



1 m p m„ 



1 + 



2m,, 



1/3 



(34) 



(35) 



Depending on the proton injected energy yo with respect to y PlP air 
and y PiP ion, the protons can produce either photopairs or both 
photopairs and photopions. We emphasise, however, that the 
above characteristic proton energies are only indicative. Since 
we are using the full emissivities, rather than delta-function ap- 
proximations for the energies of secondary particles and pho- 
tons, it is possible to have photopairs even if Eq. (f3TT > does not 
strictly hold, as relativistic protons can always pair produce with 
photons at the tail of the synchrotron distribution. 

Figure[5]shows the MW photon spectrum in the case of mo- 
noenergetic proton injection with yo - 3 x 10 5 , 3 x 10 6 , and 
3 x 10 7 , while all the other parameters have been kept con- 
stant to the values given in the previous paragraph. The run for 
yo = 3 x 10 5 marginally satisfies the photopair condition and 
not the photopion one. Thus the photon spectrum consists of the 
proton synchrotron peak and the photopair synchrotron emis- 
sion. The run for yo = 3 x 10 6 satisfies both conditions, and 
the photon spectrum shows all four features, as discussed in the 
previous section. Finally, the run for yo = 3 x 10 7 has a much 
stronger photopion than photopair component, and as a result the 
synchrotron signature of the latter lies below that of the former, 
resulting in only one peak for the synchrotron spectrum, i.e. the 
photopion one. 
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Fig. 5. Steady state MW spectrum of photons resulting from a monoen- 
ergetic proton injection with 70 = 3xl0 5 (dotted line), 3xl0 6 (full line), 



and 3 x 10 7 (dashed line). The other parameters are R 
B = 1 G, t v = 0.4, and f p csc = t a . 



3 x 10 lb cm, 



To systematise the above results, we have plotted the injected 
compactness of each individual proton process (as defined in 
Eq. (1261 1. but for (2p nj (y) corresponding to each particular pro- 
cess, rather than to total proton injection) in Fig. [6] as a function 
of jo, while all other parameters have been kept constant. For 
the injected photopion compactness we considered the sum of 
y-rays and electrons resulting from neutral and charged pion 
decay, respectively, since both contribute to the photon spec- 
trum. We note that, while the proton synchrotron term varies lin- 
early with the injected energy, the photopair and photopion terms 
show more complex behaviour. Both increase initially reaching 
a maximum, and subsequently decrease. As discussed above, 
photopair is more important for lower injected proton energies, 
but as the injection energy increases, photopion becomes dom- 
inant - for the assumed magnetic field, y p , P air = 5 x 10 5 , while 
Tp.pion = 5 x 10 6 . We also note that proton synchrotron is more 
important as a loss/injection process than the other two for all 
energies, but this depends on the value of the injected proton 
compactness. As we show in the next section this trend changes 
as the proton compactness increases. 

4.3. Increasing the injected proton compactness: From linear 
to non-linear proton cooling 

We turn next to investigating the effects that the injected proton 
compactness has on the photon spectra. In the case of proton- 
only injection like the one we are considering here, there are, up 
to a degree, profound analogies to the synchrotron - SSC rela- 
tionship of a leptonic system. There the electrons radiate syn- 
chrotron photons, and consequently upscatter them through in- 
verse Compton scattering interactions. As long as the magnetic 
energy density dominates the synchrotron photon density, the 
system can be considered in the linear regime. This situation 
changes when the synchrotron photon density dominates and the 
system becomes non-linear leading to the well-known Compton 
catastrophe. 

One can consider something analogous for the sys- 
tem we are considering at pre sent - see also §2 in 
iPetropoulou & Mastichi adisl (1201 2bl) . In hadronic systems and 
for the example we have examined above one can argue that the 
system is also in the linear regime, because protons radiate by 



py 



BH 



5.5 



6.5 

logyo 



7.5 



Fig. 6. The compactnesses of proton processes for a (5-function proton 
injection as a function of the energy y . Full line is for proton syn- 
chrotron, dashed line for photopair, and dotted line for photopion. For 
this process the sum of photon and electron injection resulting from 
neutral and charged pion decay is plotted. The other parameters are 
R = 3 x 10 16 cm, B = 1 G, l„ = 0.4, and f pesc = f cr . 



synchrotron and the thus-radiated photons are used as targets 
for photopair and photopion production. As can be seen from 
Fig- 151 the synchrotron luminosity dominates, which means that 
the cooling, however small, is regulated by this process. Figure|6] 
can also be viewed in terms of the energy-loss time scales for 
each process. Although the compactnesses depend on the proton 
density, the system's very low efficiency (as discussed in the next 
section) means that, for a given proton injection, that density will 
be practically constant. Then, for each proton injection 



ep J (r P ) = -^, 



'p,esc 



the loss terms will be 



W = T 



(36) 



(37) 



where i can be 'syn', 'BH' or photopion ('py'). Substituting in 
Eqs. (|26j! and (|29j, we find 



p,esc 



(38) 



Therefore the question that becomes relevant is what hap- 
pens to the system if the proton injection luminosity is increased 
further while the magnetic field value is kept constant. This 
would essentially mean that the photon density of the system in- 
creases over the magnetic one, and as a result the photopair and 
photopion losses/injection increase more than the proton syn- 
chrotron ones. This occurs because, while the synchrotron lumi- 
nosity depends on the proton density, the photopair and photo- 
pion luminosities depend on both the proton and photon density. 
Since the photon density depends on the proton density, we con- 
clude that the above processes depend quadratically on the pro- 
ton density provided that both use the proton synchrotron pho- 
tons as targets, which holds whenever yo > yp.pion (c.f. Eq.[35l>. 

The above can be seen in Fig. [7] which depicts the steady 
state MW spectrum in the case where the injection proton 
compactness takes the values { p -QA, 1.3, 4, 13, and 40 (bot- 
tom to top). One notices that, as the injection compactness 
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Fig. 7. Steady state MW photon spectra for (5-function proton injec- 
tion at energy y = 2.5 x 10 6 and different injection compactnesses 
i v =0.4, 1.3, 4, 13, and 40 (bottom to top). The other parameters are 
R = 3x 10 16 cm, B = 1 G, t ncsc = t cr . 



increases, the synchrotron component increases linearly while 
the photopair and photopion increase quadratically. However, 
for £ p ^40, the system undergoes a transition, and the pho- 
ton luminosity goes up by a factor of 10 4 . These types 
of transitions are a sign that the system becomes supercrit- 
ical and are caused by various feedback mechanisms - like 
the p air production - synchrotron loop dKirk & Mastichiadis 



1992) and the automatic photon qu enching one dStawarz & Kirk 



2007; Petropoulou & Mastichiadis fcOl ll) - which cause a very 
fast, non-linear proton cooling with simultaneous exponen- 
tial increase in the secondaries. The system then can either 
reach a stea dy state, as in the example shown above, or show 
limit cycles dStern & Svenssonlll99ll; iMastichiadis et al.l 120051 
Petropoul ou & Mastichiadisl 1201 2bl) . A systematic search of 
their effects are beyond the scope of the present work, but we 
will investigate their properties fully i n a forthcoming paper; see, 
however, Dimitrakoudis et alj d2012l) for some preliminary re- 
sults. For the present, it suffices to say that once the behaviour 
shown in Fig. [7] is typical, i.e. for fixed initial parameters, there 
is always a value of £ p above which the system becomes super- 
critical. We define this value as £ p , cr (R, B, y max ) and we refer to it 
in §6. A final remark concerning the above figure is that the ratio 
between the photopair and photopion injection rates, which has 
remained constant as long as the system was in the subcritical 
regime, starts changing as the system becones supercritical. This 
is because the spectrum of the target photons changes causing 
the respective rates to change accordingly. 

Figure[8]shows the injected luminosity from each proton pro- 
cess. We note that synchrotron remains linear throughout, while 
photopair and photopion are initially quadratic before they be- 
come highly non-linear once the proton density enters the super- 
critical regime. 

An interesting result of the above analysis is the following. 
In cases where y P;Pa ir < To < y P , P km, then, according to the dis- 
cussion in the previous subsection, protons can directly produce 
photopairs on proton synchrotron photons but not photopions. 
The latter, however, can be still created on the synchrotron ra- 
diation of the produced photopairs. Thus, any variation in the 
injection rate will cause the proton synchrotron photons to vary 
linearly, the radiation from photopairs quadratically (as already 
explained above), and the radiation fron photopions cubically, 
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Fig. 8. The compactnesses of proton processes for a (5-function proton 
injection as a function of the proton compactness £ p . Full line is for 
proton synchrotron, dashed line for photopair, dotted line for photopion 
- for this last process the sum of photon and electron injection resulting 
from neutral and charged pion decay is plotted. The other parameters 
are as in the previous figure. 



since their rate depends on protons and photopairs. We return to 
this point in §6. 



5. Power-law proton injection 

We next examine the case where the injected protons form a 
power law, i.e. the injection function is 



y P ,min)#(y P , 



p,max 



y) 



(39) 



where y p> mi n an d 7 P ,max are the lower and upper cut-offs of the 
proton distribution, respectively, Q Pi o is the proton normalisa- 
tion (considered, for the time being, as independent of time), 
and H(x) is the Heaviside function. From Eq.[39]one can define 
a proton compactness in complete analogy to the monoenergetic 
case - c.f. Eq.l26l 

Figure [9] depicts the spectra of photons, neutrons, and neu- 
trinos of all flavours emerging from a source with parameters 
= 10 6 , R = 3 x 10 16 cm, B = 1 G, 



10, y P ,n 



y P ,min 

s = 2, £ p — 30, and f pesc = f cr . For comparison reasons we 
have also plotted the injected proton spectrum. The neutrino 
spectrum peaks at high energies and has a distribution that re- 
sembles that of y-rays resulting from 7r°-decay. The neutrons, 
on the other hand, are more sharply peaked with a maximum 
in their distribution that is very close to y p , max . Another inter- 
esting point is that of efficiency, i.e. of the fraction of the lu- 
minosity that goes into radiation, neutrons, and neutrinos with 
respect to the total power injected into protons, which is quite 
low, at least for the parameters of the present example. Thus 
photons take about 10~ 7 of the available luminosity in protons, 
while a comparable amount goes to neutrons and neutrinos. We 
mention that for these low compactnesses, neutrons escape the 
source practically unattenuated and deposit their energy in the 
surrounding region once they are transformed back to p rotons 
dKirk & M astichiadisl ll989t iGiovanoni & Kazanasl [T990). Their 
resulting radiation will peak at UHE energies from jt° decay and 
at VHE energies from synchrotr on cooling of produced protons 
(Mastichiadis & Pro ther oell 19901) . 

In Fig. [10] which is analogous to Fig. [5] for the monoener- 
getic injection case, we treat y p max as a free parameter, taking 
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Fig. 9. Spectra of photons (small dashed line), neutrons (long dashed 
line) and neutrinos (full line) escaping from the source for initial param- 
eters: s = 2,{ p = 30, f p , csc = t a , r p ,mto = 10, r P . m ax = 10 6 , R = 3x 10 16 
cm, and B = 1 G. For comparison we have drawn the proton injected 
spectrum with a dotted line. 
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log x 



Fig. 10. Steady-state spectra of photons resulting from a power-law 
proton injection s = 2, f p = 0.3, ? p . esc = t a , y p , mm = 10, and y p-max = 
10 5 5 - 10 8 with increments of 0.5 in logarithm (bottom to top). The 
other parameters are R = 3 x 10 16 cm and B = 1 G. 

values from 10 55 to 10 8 in increments of 0.5 in logarithm. For 
Tp.max = 10 5 5 the photopair threshold is just satisfied in colli- 
sions between proton synchrotron photons and protons. Even so, 
there is a low-lumimosity photopion component created from 
collisions between synchrotron photons from photopair elec- 
trons and protons. Generally, in close analogy to the monoener- 
getic case, we find that as 7 P>max increases, photopion dominates 
the photopair; however, at very high energies of proton injec- 
tion the spectrum starts saturating as photon-photon absorption 
becomes dominant. 

In Fig.QT]we plot the steady state spectra as a function of 
the initial proton luminosity. It is analogous to Fig. [7] and one 
notices again the two basic features discussed for that figure, i.e. 
the quadratic nature of the photopair and photopion processes as 
opposed to the linear behaviour of the proton synchrotron radia- 
tion and the non-linear transition of the system to high luminosi- 
ties once the protons have reached a certain critical density. 

Figure[T2]treats the proton injection slope s as a free parame- 
ter. Since in each case the injected proton luminosity is the same, 




-10 -5 5 10 

log x 



Fig. 11. MW steady-state spectra of photons resulting from a power- 
law proton injection with y p , min = 10, y p-max = 10 6 , s = 2, r pesc = t„, 
and ( v taking the values of 0.3 to 100 with logarithmic increments of 0.5 
(bottom to top). The other parameters are R = 3 x 10 16 cm and B = 1 
G. 
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Fig. 12. MW steady-state spectra of photons resulting from variable 
power-law proton injection with y p-m i„ = 10, y p , max = 10 8 , f p = 0.1, 
fp.esc = fcr, R = 3 x 10 16 cm, B = 1 G, and s taking the values of 1.5 to 
2.5 with increments of 0.25 (top to bottom). 

injection of harder power laws mean that more of the luminos- 
ity is concentrated in higher energies, resulting in flatter spectra 
for the proton synchrotron component and also in an increase 
in the luminosity going to the photopair and photopion compo- 
nents. The latter can be seen more clearly in Fig. [13] where the 
resulting luminosities are plotted against s. Because we adopted 
a high value of the upper proton cutoff (y p , ma x = 10 s ), photopion 
dominates photopair production for flat injection spectra, while 
the two processes become comparable for steeper spectra. 

The neutrinos emitted from those same proton power laws 
at steady state are plotted in Fig. [14] Close to their upper cutoff 
they resemble power laws, which are harder than the proton ones 
by a factor of about (s - 0.5)/2.5. 

6. Time variability 

The numerical scheme introduced through Eqs. (0 to (O can 
simulate time variability by treating one or more of the free pa- 
rameters as time-dependent. Here we adopt the standard proce- 
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Fig. 13. Steady-state luminosity of photons produced by proton syn- 
chrotron (full line), electron synchrotron from pair production (long 
dashed lines), and electron synchrotron from photopion interactions 
(short dashed line), plotted against s. All other initial parameters are 
as in the previous figure. 




log (E v /m e cT 

Fig. 14. MW steady-state spectra of (electron and muon) neutrinos re- 
sulting from variable power-law proton injection as in Fig.l 121 



dure followed in such cases; i.e., we use initially constant injec- 
tion parameters until the system reaches a steady state and then 
we introduce a perturbation in one of the parameters and record 
the produced photon spectrum at each instant. 

As an indicative example, we chose first as initial parame- 
ters those that led to the steady state obtained in Fig. [4]- see also 
Fig. |7] bottom curve. We chose the particular example because, 
for this case, each hadronic component leaves a clear signature 
on the MW spectrum, and thus it will be straightforward to ex- 
amine the time response of each component. We also introduce a 
time variation to the proton injection parameter Q™ 1 in the form 
of a Lorentzian profile 



f h (t;t ,w,n) =!+(«-!) 



4(t - t ) 2 + w 2 



(40) 



The above quantity reaches a maximum for t = to, i.e. 
/l(?o; to, w, ri) = n, assuming n > 1. The quantity w is the full 
width at half maximum since /l(?o ± w;to,w, ri) — (n + l)/2. 
For the present application the parameters used were n - 3.16, 
to = 100f cr , and w = 10f cl . 



Fig. 15. Energy centered lightcurves for photon energies X\ = 1(T 3 " 85 
(dashed line), x 2 = 10~ L25 (dot-dashed line), and x 3 = 10 3 15 (short 
dashed line) resulting from a Lorentzian variation of the proton injec- 
tion parameter Q™ 1 with n = 3.2 and w = 10f cl , depicted here with full 
line. This variation was introduced on a pre-existing steady-state ob- 
tained with the parameters of figure [4] (or, equivalently, fig. [7] bottom 
curve.) For this example lines x 2 and x 3 practically coincide. 



Figure [T31 shows the resulting energy-centred (i.e. with the 
horizontal axis centred on the time when the proton energy in- 
jection becomes maximum) photon lightcurves obtained for en- 
ergies (expressed in logarithm) x\ = -3.84, X2 = -1.24, and 
X3 = 3.15. According to the discussion made following Fig. [4] 
for the initial conditions used, proton synchrotron contributes 
mainly to x\, while synchrotron from photopairs and electrons 
from photopions to x%, and electron synchrotron to x$. It be- 
comes apparent that the lightcurve at energy x\, to which proton 
synchrotron is contributing, follows the variation in the proton 
injection very closely. The other energies have a clear quadratic 
dependence on proton variation, in accordance to the discussion 
given in §4. 1 . 

Figure [16] repeats the calculations, but now we take a pro- 
ton energy injection value that satisfies the relation 7 P)P air < Jo < 
Tp.pion- The bottom full line curve that corresponds to proton syn- 
chrotron varies linearly, while the middle and top ones, which 
are at the photopair and photopion electron characteristic syn- 
chrotron energies, respectively, vary quadratically and (almost) 
cubically. The reason behind this has already been explained in 
§4.1. 

The above results are relevant when the system is in the 
subcritical regime. The situation changes drastically when the 
system is entering supercriticality. For example, Fig. [T71 shows 
the energy-centred lightcurves obtained from such a case. 
Comparing the produced lightcurves with those of the two pre- 
vious figures, we find some interesting differences, (i) Even 
though the proton compactness has been increased by a factor 
of 3, the photon luminosity increases by a factor of ~100 in all 
bands, demonstrating the high non-linearity of the system, (ii) 
The lightcurves peak about 10f cr after the peak of the perturba- 
tion, while in the previous examples the time lags were shorter, 
(iii) The relation between the three components discussed above 
is destroyed by the onset of the intense electromagnetic cascade 
caused by the supercriticality. 

To further show the complex behaviour of the system when 
it becomes supercritical, we plot in Fig. [18] a flare caused by the 
sum of two identical Lorentzians which are with peaks shifted 
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Fig. 16. A flare produced in the case with y p , pa ir < 7 P .max < Tp.pion- Full 
lines are at the characteristic energies for proton (bottom), photopair 
(middle) and photopion (top) synchrotron emission. 

by 30f cr . The ensuing flare is asymmetrical with the second one 
producing much higher luminosity than the first one. 

Therefore, it becomes evident from the above that a hadronic 
system exhibits far more complex behaviour than a leptonic one. 
The exact relationship between its main radiating components 
(proton synchrotron, photopair, and photopion) depends on both 
the injected proton energy and the proton luminosity and needs 
further investigation. Nevertheless, we present below, for the 
sake of completeness, a summary of our previous results: 

1 . Subcritical regime: For proton injection t v < { p , CI (R, B, y max ) 
the system behaves linearly; i.e., variations in proton injec- 
tion cause linear variations in the proton synchrotron com- 
ponent. The photopair component, however, always varies 
quadratically, while the photopion one varies either cubically 
(for 7 Pjlnax < r P ,pion) or quadratically (y p>max > 7 P , P ion)- In ei- 
ther case, the efficiency of the system is rather low. 

2. Supercritical regime: For ( p > ( p , cr (R, B, y max ) the system 
undergoes an abrupt phase transition. At the onset of it, even 
small perturbations in the proton injection can cause high 
amplitude variations in the secondaries. As in all dynami- 
cal systems, its exact behaviour depends sensitively and in a 
non-linear way on the initial conditions and the parameters 
of the perturbation. On the other hand, the relation between 
the various components lose the mostly quadratic depen- 
dence found in the subcritical regime owing to the electro- 
magnetic cascading, which contributes to all bands. Finally, 
the efficiency of the system can be quite high. 



7. Discussion/Summary 

Hadronic models have been used extensively for fitting the 
y-ray emission of blazars. According to their basic premises, 
y-rays can be produced either directly via proton synchrotron 
radiation or via the radiation of the secondaries resulting from 
photopion collisions. In the present paper which, in many re- 
spects, is an extension of earlier work by MK95 and MPK05, 
we have approached the model by writing and solving five cou- 
pled partial differential equations, one for protons and the rest 
for the four stable species which result from proton-photon in- 
teractions, namely electrons, photons, neutrons and neutrinos, 
thus our approach is self-consistent in the sense that it conserves 
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Fig. 17. A flare produced when the system is supercritical. The param- 
eters for the Lorentzian are as in Fig.[l5] The full lines are contributions 
from the various processes with the same order as in Fig.l 161 




i-t_o 

Fig. 18. A flare produced when the system is supercritical by the sum 
of two Lorentzians of the same amplitude. Here n = 2 while the two 
peaks occur at t\ = 85f cr and h = 115f cr respectively. For comparison 
with the previous Figure, we have chosen to = 100f cr . 

energy. However, the major improvement of the present work 
over MK95 and MPK05 was its modelling the spectra of secon- 
daries from photopion interactions (§3 ) . For this we hav e made 
systematic use of the SOPHIA results dMiicke et al.l2000l) which 
were suitably parametrised in a way so the derived distribution 
functions of secondaries can serve as source functions for their 
corresponding equations. We then modelled the proton losses 
caused by these interactions so as to maintain a detailed balance 
between the energy lost by protons and that gained by the sec- 
ondaries. This approach enables us to address topics such as the 
stability of the system, the efficiency of proton luminosity con- 
version into radiation and the relation between the underlying 
proton distribution with the emitted photon and neutrino spec- 
tra. 

Assuming a hadronic model with a monoenergetic proton in- 
jection and without any primary electron injection (§4), we have 
showed that for low proton energies the only important process, 
apart from proton synchrotron radiation, is photopair (Bethe- 
Heitler) production photopion operates, albeit at relatively low 
rates. The resulting photon spectrum shows two features, both of 
synchrotron origin, one from protons and one from the produced 
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Bethe-Heitler pairs. For even higher proton injection energies, 
photopion becomes relevant and the photon spectrum shows two 
more features, one due to the synchrotron radiation from elec- 
trons and positrons produced through charged pion decay and 
one due to y-rays from 7r°-decay (see Figs. [4] and [5J. The 
relative importance of photopair and photopion depends criti- 
cally on the injected proton energy: for low proton energies pho- 
topair dominates, while the situation is reversed for high energies 
(Fig.©. 

The above simple picture holds when the proton injected lu- 
minosity (or compactness) is low and proton cooling minimal, 
which means that the proton distribution function does not de- 
viate significantly from that at injection. In this case the pro- 
ton synchrotron luminosity is dominant while the other two pro- 
cesses contribute a smaller fraction to the total radiated lumi- 
nosity. However, as the proton compactness increases, the pho- 
topair and photopion losses increase quadratically. - note, how- 
ever, that depending on the value of the proton energy, the two 
components can vary quadratically and cubically respectively, 
while the proton synchrotron increases only linearly. This means 
that the relative contribution of the two photo-processes be- 
comes progressively more significant. Furthermore, the contri- 
bution of photon-photon cascading becomes important in redis- 
tributing the photon luminosity, thus complicating the spectral 
shape. Finally, above some critical proton luminosity, the sys- 
tem becomes supercritical and, as a result, electrons and radi- 
ation rise in an auto-regulatory way causing non-linear proton 
energy losses (see Fig.|7]and [8). 

This is a very charactersitic property of the hadronic sys- 
tems, which can lead it to high efficiencies. A detailed study 
of the system in the supercritical regime is left for a forth- 
coming publication. Preliminary calculations show that there 
are at least two loops operating that can very efficiently ex- 
tract energy from the protons and give it to leptons and pho- 
tons, namely the proton-synchrotron loop at low proton energies 
(KM92) and the automatic photon quen ching ([Stawarz & Kirkl 
120071: iPetropoulou & Mastichiadisl201 ll) at higher ones. The lat- 
ter loop was studied analytically by lPetr opoulo u & Mastichiadis] 
d2012bl) . Furthermore, it was shown that for typical blazar pa- 
rameters this loop can lead to adoption of high Doppler factors 
for the outflow 6; lower values of 6 imply higher proton lumi- 
nosities leading th e system to supercriticality and to an overpro- 
duction of X-rays (Petropoulou & Mastichiadis 2012a). 

The above results generally hold in the case where the mo- 
noenergetic proton distribution is replaced by the more astro- 
physically relevant power law (§5). In this case, the maximum 
energy of proton injection becomes an important parameter. In 
addition, results are also sensitive to the choice of the proton 
spectral slope. Interestingly enough, the calculated neutrino flux 
is flatter than the producing proton power law. This might be 
relevant to neutrino experiments, like Ice Cube. 

In the present paper we have focussed on the case where the 
soft photons required as targets for photopair and photopion in- 
teractions are provided from the synchrotron radiation of pro- 
tons themselves. While this might not be a realistic case for 
some applications - for example, in the case where there is a co- 
accelerated leptonic component, the radiated synchrotron pho- 
tons will be in direct competition with the proton synchrotron 
ones for the hadronic interactions - this approach allowed us 
to investigate the time variability in an analogous way to the 
more familiar synchrotron-SSC coupling of leptonic models. We 
found (§6) that, while the proton synchrotron part always re- 
mains linear to the proton injection, the synchrotron radiation 
from photopair, which peaks at higher energies, varies quadrat- 



ically with it. Thus the hadronic models can reproduce the 
quadratic behaviour of the SSC leptonic models and do so for the 
exact same reasons. Perhaps more interesting is that, under cer- 
tain circumstances, the photopion component that peaks at even 
higher energies (see §4) shows a cubic behaviour with respect 
to the proton injection and therefore to the proton synchrotron 
component. This is a unique feature of the hadronic models; if 
also considers the flaring behaviour exhibited when these sys- 
tems enter the supercritical regime, then it is clear that their tem- 
poral behaviour can be quite rich. 

In conclusion, we would like to emphasise two key issues 
in our approach which, in our opinion, improve the concept of 
hadronic modelling, (i) The energy conserving scheme that is in- 
troduced through the kinetic equation approach can successfully 
treat the inherent non-linearity of the system. Therefore, once 
the initial conditions are chosen, the solution of the system can 
lead either to a subcritical linear solution or to a supercritical 
nonlinear one that, occasionally, might be oscillatory. There is 
no way that the behaviour of the system can be determined a- 
priori. Therefore, time-independent approaches (e.g., adopting a 
ready distribution function for the protons) might entirely miss 
this point and lead to erroneous results, (ii) The careful mod- 
elling of photopair and photopion is crucial. We have found that 
the omission of photopair production can also lead to erroneous 
results - even in cases where the choice of initial parameters do 
not seem to favour it - when the system settles in a state where 
the internally produced soft photons from photopair production 
are non-negligible. Therefore, models that do not include it suf- 
fer a severe drawback. 

Overall, the code presented here can be a useful tool in ex- 
amining the one-zone hadronic models of active galactic nuclei 
and gamma ray bursts - the latter only in the case where the 
magnetic field is assumed to be low enough not to cause substan- 
tial synchrotron losses to pions and muons before they decay. It 
can simultaneously calculate the photon and the neutrino fluxes 
from sources with the only assumption that relativistic protons 
are injected there. It takes both photopair and photopion produc- 
tion into account with high accuracy, which are very difficult to 
model. Finally the energy-conserving scheme of our approach 
ensures that we can consider the hadronic models as dynamical 
systems and study their temporal behaviour. 
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